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Abstract. A new way to construct a fitting formula for the halo mass function is pre- 
sented. Our formula is expressed as a solution to the modified Jedamzik matrix equa- 
tion that automatically satisfies the normalization constraint. The characteristic pa- 
rameters expressed in terms of the linear shear eigenvalues are empirically determined 
by fitting the analytic formula to the numerical results from the high-resolution N-body 
simulation and found to be independent of scale, redshift and background cosmology. 
Our fitting formula with the best-fit parameters is shown to work excellently in the 
wide mass-range at various redshifts: The ratio of the analytic formula to the N-body 
results departs from unity by up to 10% and 5% over 10 11 < M/(/i _1 M Q ) < 5 x 10 15 
at z = 0, 0.5 and 1 for the FoF-halo and SO-halo cases, respectively. 

Keywords: Cosmology, Halo Mass Function 



Contents 



1 Introduction 



1 



2 Analytic framework: a review 



2 



3 Extension of the Zel'dovich model 



4 



4 Numerical tests 

4.1 Determination of the barrier heights 

4.2 Direct comparison with the N-body results 

4.3 Dependence on the filter shape 



6 

6 
10 
12 



5 Discussion 



14 



1 Introduction 

The era of precision cosmology that reigned supreme in the last century has waned, 
rendering its glory to the advent of a new era of accurate cosmology [1]. While the 
goal of the previous era was the precise measurements of a model-dependent finite 
set of the cosmological parameters, the accurate cosmology aims at addressing more 
fundamental issues such as what the origin of cosmic inflation is, what the limitation 
of the general relativity is, what drove the cosmic acceleration and etc, which requires 
to explore almost infinite parameter space spanned by all possible cosmological and 
physical scenarios. 

In precision era, resorting solely to the high-resolution simulations suffices to make 
precise theoretical predictions for the determination of the cosmological parameters. 
However, in accuracy era, given that it is extremely expensive and inefficient to sweep 
through the infinite parameter space using the simulations alone, developing robust 
analytical guidelines is highly desired to complement the numerical experiments. The 
mass function of bound halos has been spotlighted as one of those few cosmological 
probes which allow an analytical approach. 

It was ref. [2, hereafter PS] who developed for the first time an analytic formalism 
for the halo mass function under several simplified assumptions on the halo formation 
process. Later, ref. [3] introduced the excursion set theory to analytically derive the 
mass function from physical principles and found that the excursion set theory yields 
the same PS mass function for the special case of the sharp-/c space filter [see also 
4]. In the light of these landmark works, rapid progress has been made to improve 
the original excursion set formalism in two different directions. One direction was to 
refine the analytical prescriptions by incorporating into the excursion set theory more 
realistic aspects of the halo formation process such as the occurrence of the clouds-in- 
clouds, spatial correlations among the proto-halos, tidal effect from the surrounding 
matter distribution, diffusive nature of the density threshold, and etc [e.g., 5-19]. 
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The other direction focused on finding a more accurate analytic formula for the 
halo mass function by modifying the functional form of the excursion set mass function 
[20-27]. The modified functional forms are characterized by free parameters whose 
best-fit values have to be empirically determined by fitting the formulae to the high- 
resolution N-body results. Although the recent analytic formulae have been found 
to agree excellently with the N-body results, the best-fit values of their characteristic 
parameters turned out to likely vary with the mass scale, redshift and background 
cosmology [e.g., 26, 28]. 

In this paper, we present a new fitting formula for the halo mass function which 
is characterized by the free parameters independent of scale, redshift and background 
cosmology. The Jedamzik formalism [5] is adopted as our statistical framework within 
which the halo mass function is automatically normalized. Throughout this paper, 
we assume a flat ACDM cosmology with the WMAP7 parameters [30] (unless stated 
otherwise) and use the CAMB code [31] for the evaluation of the ACDM linear power 
spectrum. 



2 Analytic framework: a review 

Most of the previous formulae for the halo mass functions were constructed in the the 
following analytic framework [20-27]: 

dN(M.z) p rflno- 1 pr 

Here dN(M, z)/d In M is the differential number density of the bound halos in a 
logarithmic mass interval of [InM, InM + dlnM] at redshift z per unit volume, 
<r(M, z) = b(z)a(M, 0) is the rms fluctuation of the linear density field smoothed 
on the mass scale M at redshift z where b(z) is the linear growth factor satisfying the 
condition of b(0) — 1, p is the mean mass density of the Universe, and /(er), called the 
multiplicity function, represents a differential volume fraction occupied by those re- 
gions which satisfy a prescribed collapse condition for the halo formation. In jargon of 
the excursion set theory, a random walk proceeds as the time-like variable a increases 
until it hits a specified collapse barrier. The number of those random walks which first 
cross a given collapse barrier on scale of cx(M, z) is proportional to the number of the 
bound halos of mass M formed at z [3]. 

The functional form of f(a) that is the key quantity in equation (2.1) depends on 
the correlations among random walks as well as on the shape of the collapse barrier. 
The latter is determined by the underlying dynamics while the former is related to the 
shape of the filter used to smooth the linear density field. As mentioned in section 1, 
ref. [3] showed that the original PS mass function can be derived from the excursion 
set formalism under the assumptions that the random walks are Markovian (i.e., un- 
correlated) corresponding to the sharp-/c space filter and that the collapse barrier has 
a flat shape, 5 = 5 SC , corresponding to the spherical dynamics. The PS multiplicity 
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function derived from the excursion set theory is written as 



f PS [a(M,z)} 




71 cr(M, z) 



exp 



2cx 2 (M, z) 



(2.2) 



where the height of the spherical collapse barrier 5 SC ~ 1.686 depends only very weakly 
on the background cosmology [32]. 

As also mentioned in section 1, much effort has been made to find a better fitting 
formula by phenomenologically modifying the above PS multiplicity function. For 
instance, ref. [25] modified the PS multiplicity function into 



where the four coefficients, A, a, b and c, are free parameters whose best-fit values were 
determined empirically by comparing equation (2.3) with the N-body results. Their 
formula were found to agree with the high-resolution N-body results with accuracy up 
to maximum 10% error for the case that a bound halo is identified by the friends-of- 
friends (FoF) algorithm [33]. 

For the case that a bound halo is identified by the spherical overdensity algorithm 
[34], the following formula for the multiplicity function provided by ref. [27] has been 
found to work better [see also 23, 24]: 



The best-fit values of the four coefficients A, B, C and D were determined by comparing 
equation (2.4) to the results from the high-resolution N-body simulations for a ACDM 
cosmology with the WMAP5 parameters [35]. Their formula was demonstrated to reach 
the accuracy level up to maximum 5% error over a wide mass range of 2.4 x 10 10 < 
M / (h~ l M & ) < 10 15 . For the other fitting formulae formulated within the framework 
of equation (2.1), see Table 7 in ref. [27]. 

It is worth noting that the multiplicity function in equation (2.1) does not auto- 
matically satisfy the normalization constraint of J °° daf(a) = 1 (under the assumption 
that all initial regions would eventually collapse to form bound halos) but its overall 
amplitude had to be treated as an additional fitting parameter. For instance, the pa- 
rameter A in equation (2.3) and D in equation (2.4) had to be introduced to satisfy 
the normalization constraint. 

It has been conventionally thought that the normalization of the halo mass func- 
tion is related to the issue of taking into account the occurrence of the clouds-in-clouds 
(underdense regions embedded in larger overdense regions) and that for the special 
case of Markovian random walks with flat spherical collapse barrier the overall nor- 
malization factor of 2 introduced by PS can be justified by the excursion set theory 
[3-5] as a solution to the clouds-in-clouds problem. Very recently, however, ref. [17] 
argued that the normalization of the halo mass function in fact reflects the limitation 




(2.3) 




(2.4) 
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of the excursion set approach itself having strong dependence on the correlations of 
random walks. 

Given the downside of equation (2.1) that it requires an additional free param- 
eter to satisfy the normalization of the halo mass function, here we instead consider 
the Jedamzik equation [5] as our analytic framework which automatically yields the 
normalized mass function: 



where F(M) is the cumulative volume fraction occupied by those random walks which 
exceed the collapse barrier on the mass scale M, while P{M, M'; M'bing isolated) is 
the conditional probability that a random walk exceeds the collapse barrier on the 
mass scale M provided that it for the first time just touched the collapse barrier on 
some larger mass scale M'. That is, the isolated M' means that the random walk has 
never touched (nor exceeded) the collapse barrier on mass scales larger than M'. 

For the special case of the sharp- k space filter (Markovian random walks), the 
conditional probability P(M, M'\ M'bing isolated) can be simply given as the condi- 
tional probability P(M, M') that a random walk exceeds the collapse barrier on the 
mass scale M provided it just touched the collapse barrier on some larger mass scale 
M' [5]. For non-Markovian case (corresponding to the cases of the top-hat or Gaus- 
sian filters), however, we have inequality P(M, M'\ M' bing isolated) ^ P{M, M'), and 
thus it is expected that P{M, M'; M' bing isolated) would have much more complicated 
expression. Nontheless, in the current work, we use P{M, M') for the conditional prob- 
ability in equation (2.5), since our goal here is not to estalish a physical model but just 
to find a new fitting formula for the halo mass function which is otherwise devoid of 
physical content. 

3 Extension of the Zel'dovich model 

The original Jedamzik formalism employed the spherical collapse dynamics for which 
the collapse barrier is flat (i.e., scale-independent), expressed in terms of the linear 
density contrast, just as in the original PS formalism. It has long been realized, 
however, that the condition for a given Lagrangian region to form a bound halo should 
depend not only on its initial spherically averaged overdensity but also on the tidal 
shears from the surrounding matter distribution [e.g., 7-9, 11, 20]. The tidal shears 
have an effect of disturbing the gravitational collapse and deviating the collapse process 
from spherical symmetry. To take into account the tidal shear effect, the gravitational 
collapse process should be described by more realistic ellipsoidal dynamics rather than 
by the simplified spherical dynamics. Unlike the spherical case, however, there is no 
unique collapse condition for the case of the ellipsoidal collapse process [36]. Among 
several different ellipsoidal dynamics that were used in the literature, the Zel'dovich 
approximation [29] is one of those models which can be relatively easily implemented 
into the mass function formalism since its ellipsoidal collapse condition is simply given 




(2.5) 
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in terms of the three eigenvalues, Ai, A2, A3 (in a decreasing order, Ai > A2 > A3) of 
the initial shear tensors (T^). 

In several literatures, the Zel'dovich approximation has already been used as an 
underlying dynamics for the halo mass function. For example, refs. [8, 9] assumed that 
the gravitational collapse occurs when the largest shear eigenvalue, Ai, reaches some 
critical value, A lc , which corresponds to the collapse along the first principal axes of 
the initial shear tensors in the Zel'dovich approximation. Ref. [20] derived the halo 
mass function with ellipsoidal collapse barrier of A3 = A3 C that corresponds to the 
third principal axis collapse [see also 7] according to the dynamical guideline of the 
Zel'dovich approximation. 

Those preivous works which use the Zel'dovich approximation as a dynamical 
guideline assume that the gravitational collapse occurs when the shear eigenvalues are 
all positive. Unfortunately, however, ref. [37] have shown that one of the initial shear 
eigenvalues at a proto-halo region has different sign from the other two, which basically 
invalidates this fundamental assumption. Henceforth, those previous models based on 
the Zel'dovich approximation are not physical ones but only mechanistic prescriptions 
to determine the fitting formula for the halo mass function. Apart from the invalid 
collapse condition, the previous works based on the Zel'dovich model also suffered from 
the normalization problem which is attributed to the fact that even for the case of the 
sharp-/c space filter P(M, M'\ M' bing isolated) 7^ P(M,M') if the collapse conditions 
are expressed in terms of the shear eigenvalues. Very recently, ref. [38] have suggested 
that the normalization of the halo mass function in the Zel'dovich model can be satisfied 
if the Jedamzik formalism is used as a statistical framework. 

Adopting a similar statistical strategy to that of ref. [38], we construct an extended 
Zel'dovich model (EZL) for the halo mass function for three different cases: For the 
one-dimensional case (ID EZL), the free parameter is given as the threshold of the 
smallest shear eigenvalue: A c = {A3 C }. For the two dimensional case (2D EZL), theere 
are two free parameters given as the thresholds of the second to the largest and the 
smallest eigenvalues: A c = {A2 C , A3 C }, while for the three dimensional case (3D EZL), 
there are three free parameters given as the thresholds of all three shear eigenvalues: 
A c = {Aic, A 2c , A 3c }. We emphasize here that the thresholds of the shear eigenvalues 
are not true physical collapse barriers but just free parameters that characterize the 
halo mass function in the mechanistic Zel'dovich model. 

Expressing the free parameters in terms of the shear eigenvalues, we modify the 
Jedamzik framework as 

r°° M' dN 

F(M,Ac)=/ dM'—^—P(M,M',X c ). (3.1) 
J M p dM> 

Here F(M, A c ) is the cumulative probability that the initial shear eigenvalues exceed 
the threshold values of A c on the mass scale of M, while P(M, M', A c ) is the conditional 
probability that the shear eigenvalues on the mass scale of M exceed the thresholds of 
A c provided that the shear eigenvalues on some larger mass scale of M' just equal the 
thresholds. 



- 5 - 



In equation (3.1) the cumulative probability, F(M,X C ), can be analytically eval- 
uated by integrating the three-point probability density distribution, p(Xi, A 2 , Mc), 
which was derived by [39]: 

F{MX)= f Wt =1 d n \ tP (\;a), (3.2) 
Jc 

where n — 3, 2, 1 for the cases of the 3D, 2D, and ID EZL, respectively. Here the 
lower-bound C represents the multi-dimensional area over which the integration is 
performed. For the ID EZL, it is {Ai > A2 > A3 > A3 C }, while for the 2D and 3D EZL, 
they are {Ai > A 2 > A 2 c, Ai > A 2 > A 3 > A 3c }, and {Ai > Ai c , Ai > A 2 > A 2c , Ai > 
A 2 > A 3 > A 3c }, respectively 

The conditional probability, P{M, M', A c ), in equation (3.1) can be also obtained 
under the assumption of Markovian random walks as 

P(M, M', X c ) = P(\ > \ c \\' = \ c ) (3.3) 

c p{\' = \ c ) 

where A and A' denote the shear eigenvalues on two different mass scales M and M', 
respectively. The six-point joint probability density distribution, p(X, A'), has been 
derived analytically by [40] and [41]. For a practical calculation, we discretize equation 
(3.1) and reexpress it as a matrix product. Then we obtain the halo mass function 
dN/ d In M as a column vector by converting the matrix equation. 

It is worth noting that our ID EZL formula is very similar to that of [20] in the 
respect that both of the approaches employed the same ellipsoidal collapse condition of 
A3 = A3 C . In the former, however, the halo mass function is automatically normalized 
while in the latter the halo mass function had to be multiplied wrongly by a constant 
of 12.5 without taking into account the scale dependence of the normalization factor. 
It is also worth mentioning the key differences between our EZL model and that of ref 
[38]. The latter work used the uppler limit on the largest shear eigenvalue, Ai < Ai c 
to describe the formation of marginally bound superclusters. In the current EZL 
model for the mass function of bound halos, there is no such upper limit on the shear 
eigenvalue. In the next section we determine the best-fit values of the free parameters 
of our EZL formulae by fitting them to the numerical results from the high-resolution 
N-body simulations. 



4 Numerical tests 

4.1 Determination of the barrier heights 

For the empirical determination of the free parameters of our EZL formula, one has to 
bear in mind that the numerical results from the N-body simulations depend on which 
halo identification algorithm is used. Among the several halo identification algorithms 
that have been suggested, it is the friends-of-friends (FoF) [33] and the spherical over- 
density (SO) [34] algorithms that are most frequently used in the literatures. The 
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Figure 1. Ratios of the ID, 2D and 3D EZL mass functions of the FoF halos (green dot 
dashed, blue dashed, and red solid lines, respectively) to the N-body results given in ref. [27] 
as a function of the FoF mass at three different redshifts. 

former defines a bound halo as the collection of closely packed particles whose separa- 
tion distances are less than some prescribed linking length, while the latter defines a 
bound halo as the region in which the spherically averaged density around its density 
peak exceeds some critical value. Given that the mass functions of FoF halos have 
turned out to be different from those of the SO halos [e.g., 25], we consider both of the 
FoF and SO cases and determine separately the free parameters of the EZL formula 
for each case. 

For the FoF halo case, we compare our EZL mass function to the numerical results 
of [27, hereafter PPH10] (see eq. [2.4]). We first conduct a x 2 -ntting of the ID EZL 
formula to the PPH10 at z = by adjusting only the value of A3 C and find the best- 
fit barrier height of the ID EZL to be {A3 C = 0.41}. The top-hat filter is used to 
calculate the rms density fluctuation a(M). Figure 1 shows the ratios of the ID EZL 
mass functions with this best-fit barrier height to the numerical results of PPH10 as 
green dot-dashed lines over the wide mass range of 10 11 < Mj (h~ 1 M l3 ) < 10 15 at three 
different redshifts (z = 0, 0.5 and 1 in the left, middle and right panels, respectively). 
In each panel, the horizontal dotted line corresponds to unity. As can be seen, the 
agreement between the ID EZL and the PPH10 is not so good for each case, becoming 
worse as z increases. 

We similarly conduct a x 2 fitting of the 2D and 3D EZL mass functions to the 
PPH10 by adjusting their free parameters, respectively. The best-fit parameters for 
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Figure 2. Scatters of the randomly generated shear eigenvalues satisfying the condition of 
A3 > 0.32 in the a(M) — Ai and a(M) — A2 plane in the top and bottom panels, respectively. 



the 2D and 3D cases are found to be {A2 C = 0.547, A3 C = 0.32} and {Ai c = 0.567, A2 C = 
0.542, A3 C = 0.32}, respectively. The best-fit 2D and 3D EZL mass functions are shown 
as blue dashed and solid red lines, respectively in Figure 1. As can be seen, the 2D 
and 3D EZL mass functions show much better agreements with the PPH10 than their 
ID counterpart at each redshift. At all redshifts the 3D EZL mass function matches 
the PPH10 best in the whole mass range except for in the range of M < 10 12 h~ 1 M Q 
at z = where the 2D EZL formula works slightly better than its 3D counterpart. 
The ratio of the 2D (3D) EZL mass function to the PPH10 departs from unity by up 
to 10% (7%) over the whole mass range at the three redshifts. 

Noting that the best-fit value of A2 C for the 3D EZL case is very close to that 
for the 2D EZL case and that the 2D EZL formula works almost as well as its 3D 
counterpart, we investigate how the 2D EZL constraints on {A2,As} bias the allowed 
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Figure 3. Scatters of the randomly generated shear eigenvalues satisfying the condition of 
{A2 > 0.547, A3 > 0.32} in the a{M) — \\ plane. The horizontal dashed line corresponds to 
the 3D EZL threshold of X lc = 0.567. 

distribution of Ai. Generating randomly a set of three shear eigenvalues on the mass 
scale of M with the help of the Monte-Carlo algorithm developed by [10], and selecting 
only those random points which satisfy A 3 > 0.32, we determine the values of A 2 and Ai. 
Then we repeat the whole process on different mass scales. For the detailed description 
of the Monte-Carlo algorithm to generate randomly the three shear eigenvalues, see 
[10]. Figure 2 shows the scatter plots of the selected random points in the a(M) — Ai 
and a(M) — X2 planes in the top and bottom panels, respectively. As can be seen, 
there are low but non-zero probabilities on all mass scales that Ai and A2 exceed the 
thresholds of 0.567 and 0.542, respectively when A3 > 0.32 is satisfied. 

We repeat the same procedure but with the two constraints of A2 > 0.547 and 
A3 > 0.32 to see how biased the values of Ai are. Figure 3 plots the result in the 
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Figure 4. Same as figure 1 but using the Tinker08 formula as the numerical results to 
compare with the EZL mass functions of the SO halos. 

a(M) — Ai plane. As can be seen, the 3D EZL constraint of 0.567 is automatically 
satisfied by the cuts on A2 and A3, which explains why the best-fit value of A2 C for the 
3D EZL case is so close to that for the 2D EZL case as well as why the 2D EZL formula 
works almost as well as the 3D EZL one. 

Now, let us turn to the SO halo case for which we use the numerical results of [25, 
hereafter, Tinker08] (see eq.[2.3]) to compare our EZL mass functions with. Repeating 
the same x 2 -fitti n g procedure, we determine the best-fit free parameters of the ID, 2D, 
and 3D EZL formula for the SO-halo case as {A3 C = 0.32}, {A 2c = 0.56, A3 C = 0.32} 
{Aic = 0.56, A 2c = 0.557, A 3c = 0.32}, respectively. Figure 4 shows the same as 
figure 1 but for the SO halo case. As can be seen,the 2D and 3D EZL mass functions 
agree excellently with the Tinker08 at all three redshifts. The 2D (3D) EZL mass 
function reaches the accuracy level of maximum 5% (3%) error at all three redshifts, 
which indicates that our EZL formula works well for the FoF and SO cases alike but 
slightly better for the latter case. It is worth emphasizing here that we use the same 
constant values of A c to plot the EZL formula at all three redshifts in figures 1-4, which 
implies that the multi- dimensional barrier heights of our EZL formaula, even though 
determined empirically at z = 0, are redshift- independent. 

4.2 Direct comparison with the N-body results 

In this section we compare the EZL formula directly with the numerical results from 
two different N-body simulations. First, we use the catalog of the FoF halos identified in 
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Figure 5. (Top panel): Comparison of the EZL mass functions with the numerical results 
from the Millennium simulations at three different redshifts. The parameters of the 2D and 
3D EZL models are set at {A 2c = 0.547, A 3c = 0.32} and {Ai c = 0.567, A 2c = 0.542, A 3c = 
0.32}, respectively. In each top panel the errors bars are obtained as the standard deviation 
scatter among eight Jack-knife resamples. (Bottom panel): Ratios of the EZL mass functions 
to the N-body results. 



the Millennium simulations [42] for which the total number of particles (N p ), individual 
particle mass, linear box size (L p ) and input cosmological parameters were given as 
N p = 10 10 , m p = 8.6 x ltfh~ l M &} L p = 0.5/i _1 Gpc, VL m = 0.25, fi A = 0.75, h = 
0.73, cr 8 = 0.9, n s = 1.0, respectively. 

Binning the logarithmic masses of the FoF halos from the Millennium halo catalog 
and counting their number densities of those FoFo halos belonging to each logarithmic 
mass-bin per unit volume, we numerically determine, dN/dlnM, at z = 0, 0.5, 1. The 
top panel of figure 5 compares the numerical mass functions of the FoF halos (dots) with 
Jack-knife errors from the Millennium simulations with the 2D EZL (blue dashed lines) 
and 3D EZL (red solid lines) mass functions respectively, while its bottom panels show 
the ratios of the EZL mass functions to the Millennium results. The Jack-knife errors 
are calculated as one standard deviation scatter among eight Jack-knife resamples [38]. 
For the plots of the EZL mass functions in each panel, we set the free parameters at 
the best-fit values empirically determined in section 4.1: {A2 C = 0.547, A3 C = 0.32} 
and {Aic = 0.567, A 2c = 0.542, A 3c = 0.32} for the 2D and 3D cases, respectively. 
As can be seen, both of the 2D and 3D EZL mass functions agree very well with the 
Millennium results at all redshifts. 

We repeat the same calculations but using the catalog of the FoF halos identified 
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Figure 6. Same as figure 5 but with the numerical results from the MICE simulations. The 
same best-fit parameters that are used in figure 5 are also used for the evaluating of the 2D 
and 3D EZL mass functions 

in the MICE simulations [26]: for which the total number of particles (N p ), individual 
particle mass, linear box size (L p ) and input cosmological parameters were given as 
N p = 1024 3 , m p = 23.42 x 10 10 h~ l M &} L p = 3.072 h^Gpc, VL m = 0.3, fi A = 0.7, h = 
0.7, as = 0.8, respectively. Figure 6 shows the same as figure 5 but with the numerical 
results from the MICE simulations. As can be seen, the EZL mass functions with 
the same barrier heights show excellent agreements with the MICE mass functions at 
all three redshifts. Given that two simulations used two different values of the power 
spectrum amplitude (erg = 0.9 and 0.8 for the Millennium and MICE simulation, 
respectively), the results shown in figures 5-6 imply that the multi-dimensional barrier 
height of the EZL mass function should be independent of the background cosmology. 

4.3 Dependence on the filter shape 

In the excursion set theory, the mass assignment to a bound halo depends on the shape 
of the filter used to calculate the rms density fluctuation a(M). For the results shown 
in figures 1-6, we have used exclusively the top-hat filter. Using a different filter is 
likely to alter the shape of the EZL mass function and thus accordingly the best-fit 
values of the multi-dimensional barrier heights. To investigate the dependence of the 
best-fit values of the free parameters on the filter shape, we reevaluate the 2D EZL 
mass function but for the cases of the Gaussian and sharp-k space filter. 

Figure 7 shows the 2D EZL mass functions at z = for three different filter 
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Figure 7. Variation of the 2D EZL mass function with the filter shape. The same values of 
the free parameters that was determined in section 4.1 by using the top-hat filter is applied 
to all three filter cases. 

cases. For these plots, we use the same values of {A 2c , A 3c } determined for the case 
of the top-hat filter in section 4.1. As can be seen, for the Gaussian filter case the 
number densities of the bound halos in the high-mass section (M > 10 14 /i _1 M Q ) are 
significantly reduced while those of the low-mass halos increase. In contrast, for the 
sharp-/c space filter case the trend becomes reversed: more high-mass halos and less 
low-mass halos than for the top-hat filter case. 

To see if the success of the EZL formula is not limited to the top-hat filter case, we 
repeat the whole fitting process described in section 4.1 to determine the barrier heights 
of the 2D EZL formula but for the Gaussian and the sharp-/c space filter cases. Figures 
8 and 9 show the 2D EZL mass functions (solid line) with the best-fit parameters at 
three different redshifts for the Gaussian and sharp-Zc space filter cases, respectively. In 
each figure, the N-body results from the Millennium and the MICE simulations are also 
shown as dots (in the top and bottom panels, respectively). The best-fit parameters 
for both of the cases are different from those for the top-hat filter case. The higher 
value of \2c and lower value of A3 C for the sharp-/c filter case, while the lower values of 
\2c and for the Gaussian filter case. Anyway, even for the Gaussian and sharp-/;; 
space filter case the 2D EZL mass functions with the constant values of the parameters 
at all redshifts are in excellent agreements with the N-body results. 
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Figure 8. Same as the top panels of figures 5-6 but for the case of the Gaussian filter. 

5 Discussion 

Despite its excellent agreements with various N-body results, it is not a physical model 
but only a phenomenological fitting formula for the halo mass function since the best- 
fit values of its parameters are not related to the physics of real halo formation process. 
A critical reader might then ask what the advantage of using such a more complicated 
matrix equation as equation (3.1) for the evaluation of the halo mass function rather 
than using those simpler formulae suggested in the previous works. The merit of our 
EZL model is that it automatically satisfies the normalization constraint, as mentioned 
in 3, providing accurate fits with fewer fitting parameters than the previous models. 
Note that although the previous formulae of ref. [25] and ref. [27] have four fitting 
parameters including the normalization factor, their accuracy levels are not higher 
than our 2D EZL model which has only two parameters. Moreover, the empirically 
determined values of the free parameters of our EZL model have turned out to be 
independent of redshift and background cosmology. In our companion paper (S. Lim 

6 J. Lee 2012 in preparation), we have found that the same best-fit values of the 
parameters work even for a cosmology with primordial non-Gaussianity. 

At any rate, to develop a physical model for the halo mass function within the 



-14- 




M [h _, M ] M [h _, M ] M [h^Mg] 

Figure 9. Same as figure 8 but for the case of the sharp-fc space filter. 

framework of the Zel'dovich approximation, it will be evidently necessary to describe 
the growth of the initial shear eigenvalues by three dimensional non-Markovian random 
walks. It has been recently proved by ref. [17] that neglecting the correlations among 
random walks in the excursion set theory can significantly flaw analytic evaluation of 
the halo mass function. As an example, ref. [17] demonstrated that in the asymptotic 
limit of completely correlated random walks the excursion set theory in fact recovers 
the original PS mass function without the normalization factor of 2. Regarding the 
normalization issue that the overall amplitude of the excursion set mass function for 
the non-Markovian case is substantially lower than the numerical results, ref. [19] has 
argued that it should imply the limitation of the excursion set approach itself: Due to 
the wrong assumption that the differential mass function dN/dM on the right hand 
side of equation (2.1) can be simply related to the statistics of randomly placed cells 
on which the calculation of the differential volume fraction f(a) in the left-hand side of 
equation (2.1) is based, the excursion set mass function for the realistic non-Markovian 
case underestimates the abundance of bound halos. 

The recent work of [18] may provide a clue to how to incorporate the non- 
Markovian nature of random walks of shear eigenvalues into the Zel'dovich model. 
Ref. [18] have derived a very simple and accurate analytic model for the distribution of 
the non-Markovian random walks that first cross a given collapse barrier of arbitrary 
shape, and claimed that their model is valid even for the case that the collapse barrier 
is not a function of linear density contrast but of other initial quantities. Our future 
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work is in the direction of finding a relation of the shear eigenvalues to the physics 
of true halo formation by extending the work of [18] to the case of three dimensional 
non-Markovian random walks of shear eigenvalues [see also 43]. 
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